﻿using Alea;
using Alea.CSharp;
using Alea.Parallel;
using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Linq;
using System.Text;
using System.Threading.Tasks;

namespace GSI.Core
{
	class IDW
	{
		//反距离加权插值算法
		public static string UseCpuIDW(List<GPoint> SamplePointList, List<GPoint> meshPointList)
		{
			string output = "";
			Stopwatch sw = new Stopwatch();
			sw.Start();
			foreach (GPoint MeshPoint in meshPointList)
			{
				//当前未知点到所有已知点的距离倒数集合
				double distance = 0;
				//距离倒数之和
				foreach (GPoint SamplePoint in SamplePointList)
				{
					//这里dis即为权重
					double dis = 1.0 / (Math.Pow(MeshPoint.X - SamplePoint.X, 2) + Math.Pow(MeshPoint.Y - SamplePoint.Y, 2));
					distance += dis;
					MeshPoint.Value += dis * SamplePoint.Value;
				}
				MeshPoint.Value /= distance;
			}
			sw.Stop();
			output += $"普通计算 反距离加权插值算法 用时 {sw.ElapsedMilliseconds} ms\n";
			return output;
		}
		public static string UseParalleIDW(List<GPoint> SamplePointList, List<GPoint> meshPointList)
		{
			string output = "";
			Stopwatch sw = new Stopwatch();
			sw.Start();
			Parallel.ForEach(meshPointList, (MeshPoint) =>
			{
				//当前未知点到所有已知点的距离倒数集合
				double distance = 0;
				//距离倒数之和
				foreach (GPoint SamplePoint in SamplePointList)
				{
					//这里dis即为权重
					double dis = 1.0 / (Math.Pow(MeshPoint.X - SamplePoint.X, 2) + Math.Pow(MeshPoint.Y - SamplePoint.Y, 2));
					distance += dis;
					MeshPoint.Value += dis * SamplePoint.Value;
				}
				MeshPoint.Value /= distance;
			});
			sw.Stop();
			output += $"并行计算 反距离加权插值算法 用时 {sw.ElapsedMilliseconds} ms\n";
			return output;
		}

		[GpuManaged]
		public static string UseGPUIDW(List<GPoint> SamplePointList, List<GPoint> meshPointList)
		{

			double[] sx = SamplePointList.Select(p=>p.X).ToArray();
			double[] sy = SamplePointList.Select(p=>p.Y).ToArray();
			double[] sz = SamplePointList.Select(p=>p.Value).ToArray();
			double[] mx = meshPointList.Select(p=>p.X).ToArray();
			double[] my = meshPointList.Select(p=>p.Y).ToArray();
			double[] mv = meshPointList.Select(p=>p.Value).ToArray();

			int count = sx.Count();
			string output = "";
			Stopwatch sw = new Stopwatch();
			sw.Start();
			var gpu = Gpu.Default;
			gpu.For(0, meshPointList.Count(), (m) => {
				double distance = 0;
				for (int n = 0; n < count; n++)
				{
					double dis = 1.0 / ((mx[m] - sx[n])*(mx[m] - sx[n]) + (my[m] - sy[n])*(my[m] - sy[n]));
					distance += dis;
					//这里dis即为权重
					mv[m] += dis * sz[n];
				}
				mv[m] = mv[m]/distance;
			});
			sw.Stop();
			Parallel.For(0, meshPointList.Count, (m) =>
			{
					meshPointList[m].Value = mv[m];
			});
			output += $"GPU并行计算 反距离加权插值算法 用时 {sw.ElapsedMilliseconds} ms\n";
			return output;
		}
	}
}
